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Abstract 

We derive a multiple integral representing the ground state density matrix 
of a segment of length m of the XXZ spin chain on L lattice sites, which 
depends on L only parametrically. This allows us to treat chains of arbitrary 
finite length. Specializing to the isotropic limit of the XXX chain we show 
for small m that the multiple integrals factorize. We conjecture that this 
property holds for arbitrary m and suggest an exponential formula for the 
density matrix which involves only a double Cauchy type integral in the 
exponent. We demonstrate the efficiency of our formula by computing the 
next-to-nearest neighbour zz-correlation function for chain lengths ranging 
from two to macroscopic numbers. 
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1 Introduction 

Formally, integrable systems at finite temperature have much in common with finite- 
length systems. In the former case the free energy in the thermodynamic limit can 
be expressed as the logarithm of the dominant eigenvalue of a quantum transfer ma- 
trix [23, 24], whereas in the latter case the logarithmic derivative of the largest eigen- 
value of the usual row-to-row transfer matrix [1] determines the ground state energy 
of the system of length L. In both cases the technique of non-linear integral equa- 
tions [20,21] can be applied to express the transfer matrix eigenvalue as an integral 
over appropriately defined auxiliary functions. The integrals can be evaluated numer- 
ically yielding high precision data for thermodynamic properties at arbitrary tempera- 
tures, or for the ground state energy at arbitrary lengths, respectively. 

Here we show for the XXZ spin- 1/2 chain that this formal similarity persists for 
a multiple integral representation of the density matrix of a chain segment which was 
first derived for the ground state of the infinitely long chain [15, 16, 19] and later gen- 
eralized to finite temperature [10, 12]. We also show examples which suggest that 
the factorization of the multiple integrals, that was proven for the ground state of the 
infinite chain at vanishing magnetic field [3,6,7] and recently observed at finite tem- 
perature and non-zero magnetic field [2], might also generally hold for the ground state 
of a finite chain. 



2 The XXZ chain and its integrable structure 

The XXZ chain is an anisotropic generalization of the Heisenberg spin chain. If the 
value of all local spins is 1/2 the model is integrable, and its Hamiltonian can be ex- 
pressed through the local action of the Pauli matrices a*, a y , <3 Z on L sites of a chain, 

H = J t (°5-i°5+ °5 ^ • A! °) '°) - *)) ■ (1) 

This Hamiltonian depends on two real parameters, the exchange coupling J and an 
anisotropy parameter A. We shall consider the critical antiferromagnetic regime J > 0, 
— 1 < A < 1. Setting = diag(e 1 *,e -1 *), 4> e [0,27t], we fix the boundary conditions 
requiring that 

eo\ e \\ _ 0k ( e A e A\ ~-i 



2 21 -® 2 2 ® , ( 2 ) 

where the , a, |3 = 1 , 2, denote the gl(2) standard basis (o x = e\ + e\, <5 y = \{e\ — e\ ) , 

o z = e\ — e\). We call $ the twist angle. 4> = corresponds to the familiar periodic 
boundary conditions. 

All results in this paper rely heavily on the fact that H can be derived from the 
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well-known trigonometric i?-matrix 
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,., sh(X) ... sh(ri) 

of the six- vertex model [1]. Associating a 2 x 2 L-matrix with elements 

L$(X)=F$(X)ej* (5) 
with every lattice site we can define the monodromy matrix of the XXZ chain, 

It provides, by construction, a representation of the Yang-Baxter algebra, 

&(X-t*)(T(X)®T( M )) = (T( M )®T(X))&(X-p)^ (7) 

where R = PR if P is the transposition of the two factors in C 2 <E> C 2 . We define the 
twisted transfer matrix ?(A,) = trT(k). Then, due to (7), the function ln(? _1 (0)?(X)) 
generates a sequence of commuting local operators. The first one is proportional to the 
Hamiltonian (1), 

H = 2Jsh(r ] )d x \n(r l (0)t(X))\ x=0 , (8) 

with twisted boundary conditions (2) if we identify A = ch(r|). The critical regime 
corresponds to purely imaginary T| = iy, y e [0,7t). Because of (8) we may solve the 
eigenvalue problem of the transfer matrix rather than dealing directly with the Hamil- 
tonian. 

3 The Bethe ansatz solution for the ground state 

The twisted transfer matrix t(k) can be diagonalized by means of the algebraic Bethe 
ansatz. Since this technique has been explained elsewhere (see e.g. [22]), we may be 
content here with a mere description of the result. Eigenstates \{k}) of t(X) axe gen- 
erated by the multiple action of the operators B(k), equation (6), on the ferromagnetic 
state |0) = (lf L , 

|{X»=B(Xi-3)...B(^-5)|0>. (9) 



4 



Here the set of so-called Bethe roots {X} = {Xj} ! J =l is not arbitrary, but must be deter- 
mined from the Bethe ansatz equations 

e- 2i *sh L (X,-^) N s h(Xj-X k + T\) 
1 + 1 r, FT i J / =0. j=l,-..,N. (10) 

The transfer matrix eigenvalues corresponding to the eigenstates (9) are given by 

A(Aj " e } = \sh(X-X j + Z) + sh^ + T!) l\ sh(X-X j + *) ■ 

In the following we shall concentrate on the ground state properties of the XXZ 
chain. This brings about severe simplifications. We will be dealing with a single 
rather special solution of the Bethe ansatz equations, and we do not have to touch 
the delicate question under which circumstances the set of states (9) is complete. Let 
us fix an even length L of the chain. Then the ground state of the Hamiltonian (1) 

L/2 

is the transfer matrix eigenstate (9) with {Xj}j =l the unique real solution of (10) for 
N = L/2. The real and mutually distinct Bethe roots in this solution uniquely determine 
a meromorphic auxiliary function 

a(X) " l\ sHx-h-v) (12) 

in the complex plane. In terms of this function the ground state eigenvalue Ao(A,) 
becomes 

A (X)=(l + ^ + 3))^g|g=^|i. (13, 

It follows from the Bethe ansatz equations that Aq(X) is regular at the points Xj — ^, 
j = l,...,L/2. Extensive numerical studies moreover support the conjecture that 
Ao(A,) is non-zero inside a strip —\r\\ < ImX < 0. Our following treatment of the 
ground state of the finite-size system is based on this conjecture. It has strong imme- 
diate consequences. It implies that the function 1 + a(X) is analytic inside the strip 
— ^ < ImX < ^ and that its only zeros in this strip are the Bethe roots. This together 
with the obvious analytic and asymptotic properties of a(X) is enough to set up a set 
of functional equations for the second logarithmic derivatives of a(X) and 1 + a(X) 
which together with their known asymptotics determine a(X) uniquely [21]. Once the 
functional equations are formulated it is easy to transform them into non-linear inte- 
gral equations. This procedure is, however, non-unique. We usually work with two 
alternative forms of the non-linear integral equations. One (we call it the 'a-form') is 
convenient for theoretical purposes. We shall use it in order to write the density matrix 
as a multiple integral. The other form (we call it the 'bb-form') is useful for numerical 
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Figure 1: The canonical contour C surrounds the real axis in counterclockwise manner 
inside the strip — ^ < ImX < 



calculations. The bb-form was first derived in [21]. We show it later when we demon- 
strate the numerical efficiency of our formulae. Here is the a-form of the non-linear 
integral equation, 



Ina(X) = -2r3> + Lr|+Lln 



sh(X- 



^(?i-CQ)ln(l + a(cQ)). 
e 27i 



(14) 



The integration contour C is shown in figure 1 . It surrounds the real axis in counter- 
clockwise manner. The kernel (k) is defined as 



sh(2r|) 



ish(A,-r|) sh(X + r|) ' 



(15) 



We claim that inside the strip — < ImX < ^ the function a(k) as defined in (14) is 
the same as the function a(k) as defined in (12) and (10) with N = L/2. The ground 
state eigenvalue Aq(X) can also be expressed as an integral over a(k), 



lnA (X)=i^+ L(l7t o ^ + / ^^(X-to)ln(l + a(co)). 
2 Je 2n ? 



(16) 



This determines the ground state energy and the eigenvalues of the higher conserved 
quantities as a function of L. 

Another implication of the fact that the function 1 + a(k) is analytic inside the strip 
— < ImX < ^j- and that its only zeros in this strip are the Bethe roots is the formula 



L 



dro /(w) 



e2jril + a(a>) p x a'(kj) 



(17) 



which holds for any function f(k) analytic in the strip. This formula, when read from 
right to left, enables us to rewrite sums over the ground state Bethe roots as integrals 
over the 'canonical contour' C. It is one of the key tools in the derivation of the integral 
representation for the density matrix shown in the next section. 
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4 The density matrix as a multiple integral 

The density matrix is a means to describe a sub-system as part of a larger system 
in thermodynamic equilibrium in terms of the degrees of freedom of the sub-system. 
Typically the system is divided into two parts, interpreted as sub-system and environ- 
ment. Then one is usually interested in the limit when the sub-system is kept fixed and 
the size of the environment goes to infinity. Here we shall keep both parts finite and 
study the influence of the size of the environment on the sub-system which will be a 
segment consisting of the first m < L lattice sites of the XXZ chain. The environment 
will consist of the remaining sites. 
Let 

_H 

Pl = ^w (18) 
tre t 

be the statistical operator for the chain at temperature T. Then the density matrix of 
the sub-system consisting of the first m lattice sites is defined as 

D L (T) =tr m+L .. L p L . (19) 

By construction, the thermal average of every operator A acting non-trivially only on 
sites 1 to m can now be written as 

(A) T = tr L .. L Ap L = tr L .. m A L .. m tr m+ i..xpL = tr L „ m A L „ m D L (T) , (20) 

where A\__ m is the restriction of A to a chain consisting of sites 1 to m. In particular, 
every two-point function of local operators in the segment 1 to m of the XXZ chain 
can be brought into the above form. If we follow the common convention and use the 
same symbols for the local operators and for their restriction to the first m sites, 
we find the expression 

DL%£{T) = tri.. M ei*...e m fcD L (T) = (eft . . .e m fc) T (21) 

for the matrix elements of the density matrix. 

Here we are interested in the unique (normalized) ground state |¥o) of the system 
of finite even length. In the limit T — > 0+ the statistical operator converges to the 
projector |*Fo)(¥o| onto the ground state, and the formula (21) for the density matrix 
elements turns into 

D 4i~t = T l % D ^\:t:^ = <*okij£ (22) 

We shall use a trick suggested in [18] in order to express (22) entirely in terms of 
data related to the monodromy matrix T(k). Setting X = in the second equation (6), 
using Lj(0) = Roj(0) = Poj (if denotes the auxiliary space) and comparing both 
sides [13] we obtain 

efi = tJ-\0)l£(0)t-'(0). (23) 
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It follows that 

D 4l:t = W^o). ..r£(o)r»(o)|<p >. (24) 

In order to apply the techniques developed in [11] for the finite-temperature case we 
regularize the expression by introducing inhomogeneity parameters j = 1, . . . ,ra, 
in the following way. Define an 'inhomogeneous density matrix' 

(Wl^ 1 (^i - 5) - • • r^fe - n } | W > 



where |{A,}) is the (unnormalized) Bethe ansatz ground state. Then 



n CX1...CX, 

^ L (3i...p„ 



lim 

^1 v;^m _ 



"7 



■ft, &l 



(25) 



(26) 



The expression (25) is of the same form as in the finite-temperature case consid- 
ered in [10, 12] with the monodromy matrix elements obeying the same commutation 
relations (7). Moreover we have an auxiliary function a(k) which shares some of the 
features of the finite-temperature auxiliary function and satisfies, in particular, equa- 
tion (17). The inhomogeneous density matrix can therefore be represented as a multi- 
ple integral following the same lines of reasoning as in [10]. Since the calculations are 
very similar we can skip all details here and present the final result. 

Let |a + | be the number of up-spins (or ones) in the sequence of upper indices 
(ctj)J =l of the inhomogeneous density matrix element (25) and the number of 
down-spins (or twos) in the sequence of lower indices. The conservation of the z- 
component of the total spin implies that all density matrix elements with |oc + 1 + | (3 | 7^ 
m must vanish. Hence, 1 |3 | = m — |oc + 1 for the non-vanishing density matrix elements. 

Those are conveniently labeled by two finite sequences of positive integers an ^ 
(jyt)™ = | a +| + P where xj denotes the position of the (|a + | — j + 1 )th up-spin in (cx 7 )? =1 , 
and yk denotes the position of the (k — |a + |)th down-spin in ((3,-)™ =1 . Then 



n a i- 
Dl Pi. 



•OCm/t 



n 



7=|a 



yj- 



det(-G(co;,^)) 



(27) 



Ul<j<k<m M^k -^j) Sh(t0 7 - - t0 fe - Tl) ' 

where a = 1 /a and where the function G((0, ^) has to be calculated from the linear 
integral equation 

sh(ri) f J(oG((0,^) 



G(U) = 



+ 



sh(X-5)sh(X-£-Ti) y e 27i(l + a((o)) 



L 



■K^X-io). (28) 
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The contour C in (27), (28) is the same as in figure 1. Remarkably, (27) and (28) 
are of precisely the same form as in the finite-temperature case, the only difference 
being the definition of the auxiliary function a, equation (14). Thus, many results that 
were obtained for the finite-temperature case can be carried over to the finite-size case 
without further effort. 

Performing the homogeneous limit t,j — > ^ in (27) we obtain a multiple integral 
formula for the density matrix D^o 1 " o m . This limit was described elsewhere [12, 19]. 
Here we have to take into account that in the derivation of (27), (28) we have assumed 
that the ^ lie inside C. Thus, in order to calculate the homogeneous limit, we first 
have to push the contour to It turns out that the multiple integral formula for the 
homogeneous density matrix is not very efficient numerically (see [8] for the finite- 
temperature case). For this reason and for space limitations we leave the homogeneous 
limit of the multiple integral as an exercise to the reader. In the next section we shall 
rather proceed along the lines of the recent paper [2] where for the isotropic model the 
inhomogeneous formula (27) was first split into a sum over products of single integrals 
and where the homogeneous limit was performed only after that. To be more precise, 
such 'factorization' was carried out for m = 2, 3 and then a general formula inspired 
by [5] was conjectured for finite temperature but zero magnetic field. 



5 Factorization for XXX 



In the following we restrict ourselves to the isotropic limit A — > 1. In order to perform 
this limit in our formulae we have to replace T| by i£ with £ — > 0. In a similar way 
we have to rescale the spectral parameter X in (14), (28), the inhomogeneities, the 
integration variables and the functions a and G. Then 



n Li4..| l ncr^ i ) n 

•7=|a+|+l JL Z71 " + a \ w j)) k=l k=yj+l 
detG(co,-,^) 



ni<j<k<m&k-Z,j)(®j-®k-i) ' 
The rescaled auxiliary function a satisfies the nonlinear integral equation 

lna(X) = -2i4> + Lln^-^ M^Ml + a(co)) 



(29) 



\X+U h n l + (X-co) 2 



(30) 



and a = 1/a. The contour now surrounds the real axis counterclockwise slightly below 
ImX = j and slightly above ImX = —j. The rescaled G in the XXX limit is defined 
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by an integral equation on the same contour which reads 

1 f d(H G(o>,5) 



(3D 



¥( ^ 2) " / e 7t(TTaM) (co-^fco-^-i) • (32) 



y e ji(l + a(a))) l + (X-a)) 2 

For the following it is important to notice that (29) and (3 1) are exactly of the same 
form as in the finite-temperature case [2], where (for m = 2, 3) (31) was used in order 
to factorize (29). The argument did not depend on the form of a and therefore applies 
here in exactly the same way. We just have to replace the finite-temperature auxiliary 
function used in [2] by the finite-size auxiliary function (30). 

Let us review the results of [2]. It was shown that the inhomogeneous density 
matrix for m up to 3 can be expressed in terms of functions defined by single integrals. 
The most important one is 

dcoG(co,£i) 1 

k rc(l + a(a>)) (co-^ 2 )(to-^ 2 -i) 

We suggest (see the conjecture below) that this is the only transcendental function 
needed in the description of the inhomogeneous density matrix for arbitrary m and 
vanishing twist angle <I> = 0, and that in this case the length dependence of the density 
matrix enters through \|/(£i,^2) alone. If we want to consider non-zero <I> we have to 
deal with another family of functions 

which where called moments in [2]. 

For a compact notation of our final formulae it turns out to be useful to introduce 
certain combinations of the moments with rational functions. We first of all note that 
in the thermodynamic limit for zero twist angle the moments turn into polynomials in 
£, of order j — l, 



k=0 



(34) 



The 'normalized moments', 



9^)=<M$)-4>f(S)> (35) 
then vanish for 1/L,4> — > 0. We use them to define the symmetric combinations 

An(Sb--->Sn) — = , (36) 

lll<;<A;<n Sfcj 

where the shorthand notation h^j = — £,j was employed. The first moment cpi is 
exceptional among the (p 7 in that it becomes trivial even for finite length if only the 
twist angle vanishes, 

lim<pi(5)=0. (37) 
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It follows that 

lim A; (£) = , for all j E N. (38) 
Instead of \|/(^i,^2) we shall use the closely related expression 

y(^ 2 ) =[l + Gi-k) 2 M$i,k)-l (39) 
in terms of which our final formulae look neater. We also define 

y (^ 2 ) = Jim Y(^ 2 ) . (40) 

All density matrix elements for m = 1,2,3 can be written in terms of these func- 
tions. A complete list can be found in the appendix of [2] where the functions y(£,j,^k) 
and A„(^i, . . . have to be inserted according to our definitions above. To give 
examples let us only recall the expressions for the emptiness formation probabilities 
here, 

DL\\(lub)=\-^uh) + \(Mll)+Mb)) + lMb&), (41a) 

40i;i2<;i3 24^13^23 60 

_ 3 + 2^^13^3 perms _ (41b) 

120^13^23 

In the untwisted limit (38) applies and our result reduces to 

D L \\(^ 2 )J--^y (^ 2 ), (42a) 

D L }1}(^2^3) =^ + Vj'f 23 ttfa.k) + cyclic perms. (42b) 
24 24^13^23 

Note that the only effect of taking the limit 1/L — > here is that the function 
To (^1,^2) changes into its limiting form 

i lirn o Yo(^i,^) = 2[1 + -fe) 2 ]^, -fct) - 1 , (43) 

where 

As in the temperature case a similar statement holds true for all density matrix elements 
for m = 1, 2. The rational prefactors of Yo(^j, are the same as in the thermodynamic 
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limit. For this reason it was conjectured in [2] for the temperature case that the expo- 
nential formula for general m obtained in [5] holds also for finite temperature. Further 
evidence for this conjecture was supplied by the comparison of high-temperature ex- 
pansion data for the multiple integrals and for the conjectured exponential formula for 
m = 3,4. Regarding our results described above it seems likely that the scope of the 
exponential formula is even wider and that it also holds in the finite-length case under 
consideration. Let us briefly recall how it looks like. 

In order to obtain a convenient description of all density matrix elements we shall 
resort to a notation that we borrowed from [4]*". We arrange them into a column vector 
h m E ((C 2 )® 2 " 1 with coordinates labeled by +, — instead of 1, 2 according to the rule, 

7=1 



where Xj = — i£j for j = 1 , . . . , 



m. 



Conjecture. The density matrix of a finite sub-chain of length m of the XXX chain of 
finite length L is determined by the vector 

h m (X h ...,X m ) = — eP«^-M Smi Sm = Y[sj], (46) 

2 7=1 

fl »(^-W = — 4 ] jmm [1-( W -A* 2 ) 2 ] 2 (47) 

xtr Wj2) 2 )2 {7'(^;Xi,...,A 7B )(8) [T( Ml ;X l ,...,X m )®T(n 2 ;X l ,...,X m p-]} , 

through (45). By the integral over m, p. 2 it is meant to take the residues at the poles 
Xi, . . . , X m of the integrand. 

For the notation we are referring to [5]H; The vector s = (q) <S> (j) — (j) <8> ( ) is 
the spin singlet in C 2 <g> C 2 . The vector spaces in (C 2 )® 2m are numbered in the order 
1,2,. ..,n,h,n— 1,...,1. This defines s m . T~ is the projector onto the one-dimensional 
subspace of C 2 <g> C 2 spanned by s. 

In order to define the transfer matrices in the integrand in (47) we first of all intro- 
duce an L-matrix L(X) G U (sl 2 ) ®EndC 2 , 

L{1)= 2^5 (2X+1+Ea ®° a) ' (48) 



"This definition was first introduced in [3] and later modified in [4]. 

II In fact, the only difference between our formula (46), (47) and the result of [5] is in the function yo. 
In [5] a function co was used which is related to yo by 

0)(M - M) = hm 



i/l^o2(1-(X 1 -X 2 ) 2 )' 
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where the E a G sh are a basis satisfying [E a , E^] = 2i£ a P Y E' Y , where d is determined by 
the Casimir element through d 2 = (E a ) 2 + 1 and where p(X, d) satisfies the functional 
relation 

l-lX-d 

p(X,d)p(X-l,d)= 2X _ d (49) 
(for more details see [5]). Then, for integers, the 'transfer matrices' 

tr z r(A,;A,i,...,A, n ) = 

tx z L l {X-X l -\)...Ln{X-X n -\)L n {X-X n )...L l {X-X l ) (50) 

entering (47) are defined by substituting the irreducible representation of U(sl2) of 
dimension z into the definition (48) of the L-matrices. For non-integer z this can be 
analytically continued into the complex plane. 

6 A numerical case study 

Finally we would like to demonstrate that the formulae obtained above are numer- 
ically efficient, at least for small m. The examples we will be focusing on are the 
zz-correlation functions (o^Oj) and (o^o,) as functions of the chain length L in the 
untwisted case 4> = 0. For these it is sufficient to know the emptiness formation prob- 
ability for m = 2, 3, since 

(o\o z 2 ) = 4D L \\ - 1 , (a\al) = SD L \\\ - SD L \\ + 1 . (51) 

Hence, we have to perform the homogeneous limit £,j — > j in (42), yielding 

(°?o|> = " 5Yo(|, k) = k~ I) , (52a) 

=1 - * v ( * , i) - 1^(1, i) + I ¥ xy(i i) , (52b) 

where we denote derivatives with respect to the first and second argument by subscripts 
x and y, respectively. 

In order to calculate the functions \\f, \\fxx and ^ xy on a computer we switch to the 
bb-formulation [8,21] mentioned in section 3. For real x we define b(x) = a(x+ 5) 
and b(x) = a(x - \ ) . Then [21] 

lnb(x) =Lln(th(jw/2)) + j ^3C(x -y) ln(l + b(y)) 

- ^K(*-y + i-i0)ln(l+%)), (53a) 
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lnb(x) =Lln(th(jcc/2))+ f Jp3C(x -y) ln(l + b(y)) 

- £^X(*-y-i + iO)ln(l + b(y)), (53b) 

where %{x) is defined by (44). The function \\f and its derivatives can be expressed 
in terms of b and b. For this purpose we also have to adapt the form of the function 

G(k,Q. Following [8] we define g^\x) = ±G(x± 5,^). These functions satisfy the 
linear integral equations 

(+)/ s n r & ys ( t\y) nrl 

- ^ 3C(x-y + i-iO), (54a) 

J -°°2n(l + b (y)) 



n r°° dyg^ \y) 



ch(nfe-x)) 7—231(1 + b (y)) 



r dy^ +) (y) 

-/■. 2. ( i + V' M ) 3C(j - y - i+M)) - (54b) 

Using b, b and we can express the function \|/(£i,^2) as 



dx 

ch(7C(^2-Jf)) 



l + b- 1 ^) l + b _1 (x) 



(55) 



This formulation is now convenient for the numerical evaluation of the zz-correlation 
functions (52) for which we need \\f(j,j), ^xx{\, j) and V|/ x> (^, j). For the expansion 
of the kernel function we can use 

oo 

3CW=2£(-l)t(l+l)^, (56) 

where C#(x) = L^=i (— l) k+l /k* is the alternating zeta series**. Then 

r°° j 
=41n2+ / — - 



idx 



1 + b l (x) 1 + b 



=6^(3) + 



cbci7tch(7Uc) 
sh 2 (7i(x + i0)) 



.(-)', 



l+b-Hx) i+r'w 



*For Rex > 1 the Riemann zeta function C,(x) and ^ fl (x) are related by t, a {x) = (1 — 2 1 
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Figure 2: Nearest neighbour zz-correlation function for chains of even length L, solid 
curve represents the analytic continuation to arbitrary real and positive L as defined by 
our integral representation. 



V „(i,i)=-6«3) + 



id* 



sh(jt(jc + iO)) 



+ 



1 + b l + b *(*). 



(57) 



where the primes denote derivatives with respect to Using (54) and (57) in (52) we 
calculated the nearest and next-to-nearest neighbour zz-correlators numerically. The 
non-linear integral equations (53) as well as the linear integral equations (54) were 
solved iteratively in Fourier space. The derivatives of gf^ with respect to £ were com- 
puted by solving the integral equations obtained from (54) by taking the derivatives 
with respect to 

Figure 2 shows the nearest neighbour correlator and figure 3 the next-to-nearest 
neighbour correlator. Since our model represents an antiferromagnet the former must 
be negative and the latter positive. The weakening of the correlation with growing 
length can be attributed to what is called 'quantum frustration' in condensed matter 
physics. Let us illustrate this notion with an example. The (unnormalized) ground 
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Figure 3: Next-to-nearest neighbour zz-correlation function for chains of even length 
L, solid curve represents the analytic continuation to arbitrary real and positive L as 
defined by our integral representation. 

states of the XXX chain of length 2 and 4, respectively, are 

|^)2 = Ut)-ltl) , 
Neel 

\ g s) 4 =2| mt) +2\ TIT!) - 1 UTT) - 1 TUT) - 1 TTII) - 1 ITT!) • (58) 

N v ' N v ' 

Neel QM frustration 

The states with alternating up and down spins on consecutive sites are called Neel 
states. They are a 'classical caricature' of an antiferromagnet. The zz-correlators in a 
Neel state are (o|a^ +1 ) = ( — 1)" and realize 'perfect antiferromagnetic order'. Such 
type of order is realized in the ground state of the XXX chain only for L = 2, where the 
correlations look like in the classical case (see figure 2) and where the next-to-nearest 
neighbour correlator is not defined. For L = 4, as can be seen from the ground state 
wave function (58), there is a certain probability to have parallel spins on neighbouring 
sites and antiparallel spins on next-to-nearest neighbour sites. This reduces the correla- 
tions in both cases. For growing chain length the 'Neel order' is even more frustrated, 
e.g. three or four parallel spins appear, and the correlations are further reduced, which 
explains the monotonous behaviour of our curves in figures 2 and 3. 
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L 


(<5\<5 Z 2 ) 




2 


-1.00000000000000 


— 


4 


-0.66666666666667 


0.33333333333333 


8 


-0.60851556815620 


0.26103720534839 


16 


-0.59519136338473 


0.24696584167998 


32 


-0.59193864328956 


0.24374937989865 


64 


-0.59113127886152 


0.24297329183505 


128 


-0.59092994011745 


0.24278223127753 


256 


-0.59087965782193 


0.24273481483257 


512 


-0.59086709385781 


0.24272300601642 


1024 


-0.59086395383499 


0.24272006021644 


CO 


-0.59086290741326 


0.24271907982574 



Table 1 : zz-correlators as functions of the system size. 



It follows by inspection of equation (53) that b(x) and b(x) vanish in the thermo- 
dynamic limit L — > oo. Hence, the integrals in (57) all vanish and V|/(j, j) = 4 In 2, 
Vxy{\, \) = 6^(3), VxjcCj, 5) = —6^(3). Inserting this into (52) we obtain 

lim (o*of) =i-f ln2 

lim (o*o§)=i-^ln2 + 3£(3) (59) 

The first equation is a corollary to Hulthen's classical result [14] on the ground state 
energy per site of the XXX chain, and the second one is a well known result due to 
Takahashi [25] which was reproduced from the multiple integral formula of Jimbo et 
al. [15] by Boos and Korepin [7]. Here we have calculated the correlation functions for 
finite chain length as correction to the asymptotic values (59). The zz-correlators for 
L = 4, far away from these asymptotic values, can be easily obtained from the ground 
state wavefunction \gs)4 in (58), (o\o z 2 ) = ~\ an d (^of) = 5- It is remarkable that 
these values are reproduced from our integral equations to 13 digits precision without 
too much effort (see table 1). 

7 Conclusions 

We have obtained a multiple integral formula (27) for the zero temperature limit of 
the density matrix of a finite segment of the XXZ chain which holds for every even 
chain length L and is again of the same form as the formerly known formula for finite 
temperature but infinite length. The multiple integrals are 'parameterized' by a pair of 
functions a(X), G(A,,£) which fix their physical meaning. For one such pair we obtain 
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the finite-temperature density matrix for another pair the ground state density matrix 
for the finite chain, the only difference being the driving term in the non-linear integral 
equation (14) for a(k). In the thermodynamic limit for vanishing twist angle (the zero 
temperature limit for vanishing magnetic field) the integrands in the multiple integrals 
turn into explicit functions and the formulae of Jimbo et al. [15, 16] are recovered. 

Even the linear integral equations for G(k, are of the same form in the tem- 
perature case and in the finite-length case. Since only this form was relevant for the 
reduction of the multiple integrals to sums over products of single integrals (the 'fac- 
torization') in the finite-temperature case [2], a similar factorized form of the short 
range correlations is valid for finite length and the conjecture formulated for the gen- 
eral finite-temperature case in [2] is likely to hold in the finite-length case as well (see 
(46)). The phenomenon of factorization of correlation functions, first observed in [7] 
and rather well understood in the thermodynamic limit for zero temperature and van- 
ishing external field by now [5], may turn out to be valid in a much broader context 
(compare the discussion in the summary of [9]) and may even turn out to be a general 
characteristic of quantum integrable models related to the Yang-Baxter equation. 

In this work we have concentrated on the density matrix, since we wanted to test 
if the factorization scheme of [2] also works in the finite-length case. We have seen 
that this is indeed the case. Concerning the multiple integral representation we have 
no doubt that similar formulae as derived for the two-point functions [10] and for a 
generating function of the zz-correlation functions [11] for finite temperatures in the 
thermodynamic limit also hold in the finite length case. It will be interesting to com- 
pare the formulae for the generating function obtainable by using the auxiliary function 
a and our function G with the result of [17], where another multiple integral for the 
finite-length system was derived. 
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